# This file runs code that reproduces all figures in the text and
# appendix. It loads all three datasets. You will need to download
# the datasets, place them in a folder, and set that folder as
# your working directory. 

# set working directory on your machine
setwd("YOUR_FOLDER") 

# load packages
library(tidyverse)
library(datawizard)
library(Rmisc)
library(broom)
library(jtools)
library(texreg)
library(stargazer)
library(gtsummary)
library(gt)

#### IMPORTING RURAL ID/INTERNET INSECURITY DATA SET ####
d <- read_csv("2021_CES.csv") 

druralDesc <- d
d <- normalize(d, 
               select = c("similarity", "harm.scale", "pps", "attr.isp", "govt.spend", "ser.imp",
                          "widespread", "harm.economic", "harm.psych", "harm", "harmxscale",
                          "serious", "important", "concerned", "faminc_new",
                          "pew_religimp", "pid7", "ideology", "age", "pro.infra"))

drural <- d



#### IMPORTING PARTY ID/AIR POLLUTION DATA SET ####
d <- read_csv("2022_Prolific.csv")

dpidDesc <- d

d <- normalize(d, select = c("pps", "policy.filters",
                                         "similarity", "cause.internal",
                                         "cause.residents", "blame.victim",
                                         "AP", "r.concern", "age", "ideo7", "pid7",
                                         "income", "religiosity", "trustpeople",
                                         "trustgov", "proximity", "cause.govt"))

dpid <- d


#### IMPORTING CYBER-ATTACK DATASET ####
d <- read_csv("2023_CES.csv")

d <- normalize(d, select = c("similarity", "pps", "policy", "blame_victim", "r_concern"))

dcyber <- d


#### TREATMENT EFFECT SIZES ####

# proportion of rural respondents who are Dem., Ind., Rep.
prop.table(table(drural$rural, drural$pid3new), 1)
 
# policy attitudes
summary(lm(govt.spend ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural))
# effect size
0.0962887 / sd(drural$govt.spend, na.rm = T)


# effect size
0.109232 / sd(drural$concerned, na.rm = T)

# policy attitudes
summary(lm(policy.filters ~ condition_inparty_w_NP + coll.degree + religiosity + trustpeople  +
                 trustgov + nonwhite + age + female_or_nb + newsint, 
           dpid))
# effect size
-0.037063 / sd(dpid$policy.filters, na.rm = T)


# personal concern
summary(lm(r.concern ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity + trustpeople  +
                 trustgov + nonwhite + age + female_or_nb + newsint, 
           dpid))
# effect size
-0.060101 / sd(dpid$r.concern, na.rm = T)



#### MAIN FIGURES ####

##### Similarity #####

s1sim <- lm(similarity ~ rural.victims*rural + pro.infra +
                  female + educ + interest +
                  age + pew_religimp, drural)

s2sim <- lm(similarity ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity + trustpeople  +
                  trustgov + nonwhite + age + female_or_nb + newsint, 
            dpid)

s3sim <- lm(similarity ~ treat_char + female + educ + interest +
                  age + pew_religimp, dcyber)

# S1 similarity
s1sim_dat <- make_predictions(s1sim, pred = "rural", pred.values = c(0,1), 
                               at = list(rural.victims = c(0,1))) %>%
      mutate(rural.victims = ifelse(rural.victims == 0, "Urban", "Rural"))



ggplot(s1sim_dat, aes(fct_relevel(rural.victims, "Urban"), 
                       similarity, group = factor(rural))) + 
      geom_line(aes(color = factor(rural)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(rural)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Similarity")


# S2: Similarity
s2sim_dat <- make_predictions(s2sim, pred = "pid3", 
                              pred.values = c("Republican", "Democratic"), 
                              at = list(condition_inparty_w_NP = c("In-party", "0No party", "Out-party")))

ggplot(s2sim_dat, aes(fct_relevel(condition_inparty_w_NP, "In-party"), 
                      similarity, 
                      group = factor(pid3))) + 
      geom_line(aes(color = factor(pid3)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(pid3)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      scale_x_discrete(labels = c("In-Party", "No Party", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Similarity")

# S3: pps ~ effect plot
s3sim_dat <- make_predictions(s3sim, pred = "treat_char",
                              pred.values = c("Out-party" , "0Control", "Primary Identity"))

ggplot(s3sim_dat, aes(fct_relevel(treat_char, "Out-party"), 
                      similarity)) + 
      geom_line(group = 1) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax), 
                      position = position_dodge(width = .2)) +
      scale_x_discrete(labels = c("Out-party", "No ID", "Primary ID")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Similarity")



##### PPS #####
s1harm <- lm(harm.economic ~ rural.victims*rural + pro.infra +
                   female + educ + interest +
                   age + pew_religimp, drural)

s2pps <- lm(pps ~ condition_inparty_w_NP + coll.degree + religiosity + 
                   age + female_or_nb + newsint, 
            dpid)

s2ppsgop <- lm(pps ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                   age + female_or_nb + newsint, 
               dpid)

s3pps <- lm(pps ~ treat_char + female + educ + interest +
                     age + pew_religimp, dcyber)

s3ppsgop <- lm(pps ~ treat_char*republican + female + educ + interest +
                  age + pew_religimp, dcyber)



# S1: economic harm ~ T*rural interaction plot
s1harm_dat <- make_predictions(s1harm, pred = "rural", pred.values = c(0,1), 
                 at = list(rural.victims = c(0,1))) %>%
      mutate(rural.victims = ifelse(rural.victims == 0, "Urban", "Rural"))

ggplot(s1harm_dat, aes(fct_relevel(rural.victims, "Urban"), 
                       harm.economic, group = factor(rural))) + 
      geom_line(aes(color = factor(rural)), 
                    position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(rural)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Economic Harm")


# S2: pps ~ T*pid3 interaction plot
s2pps_dat <- make_predictions(s2ppsgop, pred = "pid3", 
                  pred.values = c("Republican", "Democratic"), 
                  at = list(condition_inparty_w_NP = c("In-party", "0No party", "Out-party")))

ggplot(s2pps_dat, aes(fct_relevel(condition_inparty_w_NP, "In-party"), 
                      pps, 
                      group = factor(pid3))) + 
      geom_line(aes(color = factor(pid3)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(pid3)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      scale_x_discrete(labels = c("In-Party", "No Party", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Problem Seriousness")

# S3: pps ~ effect plot
s3pps_dat <- make_predictions(s3pps, pred = "treat_char",
                              pred.values = c("Out-party" , "0Control", "Primary Identity"))

ggplot(s3pps_dat, aes(fct_relevel(treat_char, "Primary Identity"), 
                      pps)) + 
      geom_line(group = 1) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax), 
                      position = position_dodge(width = .2)) +
      scale_x_discrete(labels = c("Primary ID", "No ID", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Perceived Problem Seriousness")

# effect size in SDs of DV
-0.0464832 / sd(dcyber$pps, na.rm = T)



##### Policy Attitudes #####

s1pa <- lm(govt.spend ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural)

s2pa <- lm(policy.filters ~ condition_inparty_w_NP + coll.degree + religiosity +
                 age + female_or_nb + newsint, 
           dpid)

s2pagop <- lm(policy.filters ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity + 
                     age + female_or_nb + newsint, 
              dpid)

s3pa <- lm(policy ~ treat_char + female + educ + interest +
                  age + pew_religimp, dcyber)


# S1: pa effect plot
s1pa_dat <- make_predictions(s1pa, pred = "rural", pred.values = c(0,1), 
                                           at = list(rural.victims = c(0,1))) %>%
      mutate(rural.victims = ifelse(rural.victims == 0, "Urban", "Rural"))

ggplot(s1pa_dat, aes(fct_relevel(rural.victims, "Urban"),
                     govt.spend, group = factor(rural))) + 
      geom_line(aes(color = factor(rural)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(rural)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) + 
      labs(x = "Victims' Identity", y = "Policy Support")

# S2: pa effect plot by party
s2pa_dat <- make_predictions(s2pagop, pred = "pid3", 
                              pred.values = c("Republican", "Democratic"), 
                              at = list(condition_inparty_w_NP = c("In-party", "0No party", "Out-party")))

ggplot(s2pa_dat, aes(fct_relevel(condition_inparty_w_NP, "In-party"), 
                      policy.filters, 
                      group = factor(pid3))) + 
      geom_line(aes(color = factor(pid3)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(pid3)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      scale_x_discrete(labels = c("In-Party", "No Party", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Policy Support")

# S3: pps ~ effect plot
s3pa_dat <- make_predictions(s3pa, pred = "treat_char",
                              pred.values = c("Out-party" , "0Control", "Primary Identity"))

ggplot(s3pa_dat, aes(fct_relevel(treat_char, "Primary Identity"), 
                      policy)) + 
      geom_line(group = 1) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax), 
                      position = position_dodge(width = .2)) +
      scale_x_discrete(labels = c("Primary ID", "No ID", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Policy Support")


##### Personal Concern #####

s1pc <- lm(concerned ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural)

s2pc <- lm(r.concern ~ condition_inparty_w_NP + coll.degree + religiosity +
                  age + female_or_nb + newsint, 
           dpid)

s2pcgop <- lm(r.concern ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                    age + female_or_nb + newsint, 
              dpid)

s3pc <- lm(r_concern ~ treat_char + female + educ + interest +
                 age + pew_religimp, dcyber)

# S1: pc
s1pc_dat <- make_predictions(s1pc, pred = "rural", pred.values = c(0,1), 
                             at = list(rural.victims = c(0,1))) %>%
      mutate(rural.victims = ifelse(rural.victims == 0, "Urban", "Rural"))

ggplot(s1pc_dat, aes(fct_relevel(rural.victims, "Urban"),
                     concerned, group = factor(rural))) + 
      geom_line(aes(color = factor(rural)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(rural)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Personal Concern")

# S2: pc of GOP vs. Dems

s2pc_dat <- make_predictions(s2pcgop, pred = "pid3", 
            pred.values = c("Republican", "Democratic"), 
            at = list(condition_inparty_w_NP = c("In-party", "0No party", "Out-party")))

ggplot(s2pc_dat, aes(fct_relevel(condition_inparty_w_NP, "In-party"), 
                      r.concern, 
                      group = factor(pid3))) + 
      geom_line(aes(color = factor(pid3)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(pid3)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      scale_x_discrete(labels = c("In-Party", "No Party", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Personal Concern")


##### Causal Attribution #####
s1ca <- lm(attr.isp ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural)

s2cagop <- lm(cause.residents ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity + 
                     age + female_or_nb + newsint, 
              dpid)

s3ca <- lm(blame_victim ~ treat_char + female + educ + interest +
                 age + pew_religimp, dcyber)

# S1: ca
s1ca_dat <- make_predictions(s1ca, pred = "rural", pred.values = c(0,1), 
                             at = list(rural.victims = c(0,1))) %>%
      mutate(rural.victims = ifelse(rural.victims == 0, "Urban", "Rural"))

ggplot(s1ca_dat, aes(fct_relevel(rural.victims, "Urban"),
                     attr.isp, group = factor(rural))) + 
      geom_line(aes(color = factor(rural)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(rural)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "External Attribution")

# s2: ca
s2ca_dat <- make_predictions(s2cagop, pred = "pid3", 
                             pred.values = c("Republican", "Democratic"), 
                             at = list(condition_inparty_w_NP = c("In-party", "0No party", "Out-party")))

ggplot(s2ca_dat, aes(fct_relevel(condition_inparty_w_NP, "In-party"), 
                     cause.residents, 
                     group = factor(pid3))) + 
      geom_line(aes(color = factor(pid3)), 
                position = position_dodge(width = .2)) +
      geom_pointrange(aes(ymin = ymin, ymax = ymax, color = factor(pid3)), 
                      position = position_dodge(width = .2)) +
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      scale_x_discrete(labels = c("In-Party", "No Party", "Out-Party")) +
      theme_bw() +
      theme(legend.position = c(.9, .7),
            legend.text = element_text(size=12),
            axis.text = element_text(size=12),
            axis.title = element_text(size=12)) +
      labs(x = "Victims' Identity", y = "Internal Attribution")


#### APPENDIX TABLES ####

##### PPS #####

s1harmbase <- lm(harm.economic ~ rural.victims*rural, drural)

s1harm <- lm(harm.economic ~ rural.victims*rural + pro.infra +
                   female + educ + interest +
                   age + pew_religimp, drural)

s2ppsgopbase <- lm(pps ~ condition_inparty_w_NP*pid3, 
               dpid)

s2ppsgop <- lm(pps ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                     age + female_or_nb + newsint, 
               dpid)

s3ppsbase <- lm(pps ~ treat_char, dcyber)

s3ppsgop <- lm(pps ~ treat_char*republican, dcyber)


# s1 pps
screenreg(list(s1harmbase, s1harm),
       single.row = T,
       caption = "Study 1: Experimental Effects on Perceived Problem Seriousness",
       caption.above = T,
       custom.coef.names = c("Intercept", "Rural condition", "Rural", 
                             "Rural condition x Rural", "Pro-infrastructure", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s1_pps",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s2 pps
screenreg(list(s2ppsgopbase, s2ppsgop),
       single.row = T,
       caption = "Study 2: Experimental Effects on Perceived Problem Seriousness",
       caption.above = T,
       custom.coef.names = c("Intercept", "In-party", "Out-party", "Republican", 
                             "In-party x Republican", "Out-party x Republican",
                             "College degree", "Religiosity",
                             "Age", "Female", "Interest"),
       label = "app_tab_s2_pps",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, booktabs = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s3 pps
screenreg(list(s3ppsbase, s3pps),
       single.row = T,
       caption = "Study 3: Experimental Effects on Perceived Problem Seriousness",
       caption.above = T,
       custom.coef.names = c("Intercept", "Out-party condition", "Primary Identity Condition", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s3_pps",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)



##### POLICY ATTITUDES #####

s1pabase <- lm(govt.spend ~ rural.victims*rural, drural)

s1pa <- lm(govt.spend ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural)

s2pabase <- lm(policy.filters ~ condition_inparty_w_NP*pid3, 
               dpid)

s2pagop <- lm(policy.filters ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                    age + female_or_nb + newsint, 
              dpid)

s3pabase <- lm(policy ~ treat_char, dcyber)

s3pa <- lm(policy ~ treat_char + female + educ + interest +
                 age + pew_religimp, dcyber)

# s1 pa
screenreg(list(s1pabase, s1pa),
       single.row = T,
       caption = "Study 1: Experimental Effects on Policy Attitudes",
       caption.above = T,
       custom.coef.names = c("Intercept", "Rural condition", "Rural", 
                             "Rural condition x Rural", "Pro-infrastructure", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s1_pa",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s2 pa
screenreg(list(s2pabase, s2pagop),
       single.row = T,
       caption = "Study 2: Experimental Effects on Policy Attitudes",
       caption.above = T,
       custom.coef.names = c("Intercept", "In-party", "Out-party", "Republican", 
                             "In-party x Republican", "Out-party x Republican",
                             "College degree", "Religiosity", 
                             "Age", "Female", "Interest"),
       label = "app_tab_s2_pa",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, booktabs = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s3 pps
screenreg(list(s3pabase, s3pa),
       single.row = T,
       caption = "Study 3: Experimental Effects on Policy Attitudes",
       caption.above = T,
       custom.coef.names = c("Intercept", "Out-party condition", "Primary Identity Condition", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s3_pa",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)


##### CAUSAL ATTRIBUTION #####

s1cabase <- lm(attr.isp ~ rural.victims*rural, drural)

s1ca <- lm(attr.isp ~ rural.victims*rural + pro.infra +
                 female + educ + interest +
                 age + pew_religimp, drural)

s2cabase <- lm(cause.residents ~ condition_inparty_w_NP*pid3, 
              dpid)

s2cagop <- lm(cause.residents ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                    age + female_or_nb + newsint, 
              dpid)

s3cabase <- lm(blame_victim ~ treat_char, dcyber)

s3ca <- lm(blame_victim ~ treat_char + female + educ + interest +
                 age + pew_religimp, dcyber)
# s1 ca
screenreg(list(s1cabase, s1ca),
       single.row = T,
       caption = "Study 1: Experimental Effects on Causal Attributions",
       caption.above = T,
       custom.coef.names = c("Intercept", "Rural condition", "Rural", 
                             "Rural condition x Rural", "Pro-infrastructure", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s1_pa",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s2 ca
screenreg(list(s2cabase, s2cagop),
       single.row = T,
       caption = "Study 2: Experimental Effects on Causal Attributions",
       caption.above = T,
       custom.coef.names = c("Intercept", "In-party", "Out-party", "Republican", 
                             "In-party x Republican", "Out-party x Republican",
                             "College degree", "Religiosity", 
                              "Age", "Female", "Interest"),
       label = "app_tab_s2_pa",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, booktabs = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s3 ca
screenreg(list(s3cabase, s3ca),
       single.row = T,
       caption = "Study 3: Experimental Effects on Causal Attributions",
       caption.above = T,
       custom.coef.names = c("Intercept", "Out-party condition", "Primary Identity Condition", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s3_ca",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)


##### PERSONAL CONCERN #####

s1pcbase <- lm(concerned ~ rural.victims*rural, drural)

s1pc <- lm(concerned ~ rural.victims*rural + pro.infra +
                   female + educ + interest +
                   age + pew_religimp, drural)

s2pcgopbase <- lm(r.concern ~ condition_inparty_w_NP*pid3, 
                   dpid)

s2pcgop <- lm(r.concern ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                     age + female_or_nb + newsint, 
               dpid)

s3pcbase <- lm(r_concern ~ treat_char, dcyber)

s3pc <- lm(r_concern ~ treat_char + female + educ + interest +
                 age + pew_religimp, dcyber)

# s1 pc
screenreg(list(s1pcbase, s1pc),
       single.row = T,
       caption = "Study 1: Experimental Effects on Personal Concern",
       caption.above = T,
       custom.coef.names = c("Intercept", "Rural condition", "Rural", 
                             "Rural condition x Rural", "Pro-infrastructure", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s1_pc",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s2 pc
screenreg(list(s2pcgopbase, s2pcgop),
       single.row = T,
       caption = "Study 2: Experimental Effects on Personal Concern",
       caption.above = T,
       custom.coef.names = c("Intercept", "In-party", "Out-party", "Republican", 
                             "In-party x Republican", "Out-party x Republican",
                             "College degree", "Religiosity", 
                             "Age", "Female", "Interest"),
       label = "app_tab_s2_pc",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, booktabs = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s3 pc
screenreg(list(s3pcbase, s3pc),
       single.row = T,
       caption = "Study 3: Experimental Effects on Personal Concern",
       caption.above = T,
       custom.coef.names = c("Intercept", "Out-party condition", "Primary Identity Condition", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s3_pc",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

##### Similarity #####
s1simbase <- lm(similarity ~ rural.victims*rural, drural)

s1sim <- lm(similarity ~ rural.victims*rural + pro.infra +
                   female + educ + interest +
                   age + pew_religimp, drural)

s2simgopbase <- lm(similarity ~ condition_inparty_w_NP*pid3, 
                   dpid)

s2simgop <- lm(similarity ~ condition_inparty_w_NP*pid3 + coll.degree + religiosity +
                      age + female_or_nb + newsint, 
               dpid)

s3simbase <- lm(similarity ~ treat_char, dcyber)

s3sim <- lm(similarity ~ treat_char + female + educ + interest +
                  age + pew_religimp, dcyber)

# s1 sim
screenreg(list(s1simbase, s1sim),
       single.row = T,
       caption = "Study 1: Experimental Effects on Perceived Similarity",
       caption.above = T,
       custom.coef.names = c("Intercept", "Rural condition", "Rural", 
                             "Rural condition x Rural", "Pro-infrastructure", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s1_sim",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s2 sim
screenreg(list(s2simgopbase, s2simgop),
       single.row = T,
       caption = "Study 2: Experimental Effects on Perceived Similarity",
       caption.above = T,
       custom.coef.names = c("Intercept", "In-party", "Out-party", "Republican", 
                             "In-party x Republican", "Out-party x Republican",
                             "College degree", "Religiosity",
                              "Age", "Female", "Interest"),
       label = "app_tab_s2_sim",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, booktabs = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

# s3 pc
screenreg(list(s3simbase, s3sim),
       single.row = T,
       caption = "Study 3: Experimental Effects on Perceived Similarity",
       caption.above = T,
       custom.coef.names = c("Intercept", "Out-party condition", "Primary Identity Condition", 
                             "Female", "Education", "Interest", "Age", "Religiosity"),
       label = "app_tab_s3_sim",
       digits = 3,
       symbol = "+",
       float.pos = "h!", dcolumn = T, use.packages = F,
       stars = c(0.001, 0.01, 0.05, 0.1),
       include.rsquared = F)

###### strong primary IDs ######
t.test(similarity ~ treat, d, subset = treat > -1 & pi_strong == 1)
ttt <- tidy(t.test(similarity ~ treat, d, subset = treat > -1 & pi_strong == 1))
ttt <- dplyr::select(ttt, estimate1, estimate2, p.value)
stargazer(ttt, summary = F, notes = "p-value")

##### blame criterion validity #####
cor.test(drural$attr.isp, drural$pps)
cor.test(dpid$cause.residents, dpid$pps)
cor.test(dcyber$blame_victim, dcyber$pps)

cor.test(drural$attr.isp, drural$govt.spend)
cor.test(dpid$cause.residents, dpid$policy)
cor.test(dcyber$blame_victim, dcyber$policy)

##### S1 CHARACTERISTICS #####
tbl1 <- druralDesc |> select(female, age, educ, rural, race, pid3new, interest, pew_religimp) |>
      mutate(pid3new = case_when(
            pid3new == 1 ~ "Democratic", 
            pid3new == 2 ~ "Independent",
            pid3new == 3 ~ "Republican",
      ),
      race = case_when(
            race == 1 ~ "White",
            race == 2 ~ "Black",
            race == 3 ~ "Hispanic",
            race == 4 ~ "Asian",
            race > 4 ~ "Other"
      ),
      educ = case_when(
            educ == 1 ~ "No HS",
            educ == 2 ~ "HS degree",
            educ == 3 ~ "Some college",
            educ == 4 ~ "2-year",
            educ == 5 ~ "4-year",
            educ == 6 ~ "Post-grad")) |>
      tbl_summary(type = list(interest ~ "continuous",
                              pew_religimp ~ "continuous",
                              pid3new ~ "categorical"),
                  label = list(female ~ "Female", 
                               age ~ "Age",
                               rural ~ "Rural",
                               educ ~ "Education",
                               pid3new ~ "Party ID",
                               race ~ "Race",
                               interest ~ "Political interest (1-4)",
                               pew_religimp ~ "Religiosity (1-4)"),
            statistic = list(
            all_continuous() ~ "{mean} ({sd})",
            all_categorical() ~ "{p}%"),
            missing = "no") |>
      modify_header(label = "Variable") |>
      bold_labels()

##### S2 CHARACTERISTICS #####
tbl2 <- dpidDesc |> select(female_or_nb, age, educ, rural,
                           contains("race"), pid3, newsint, religiosity) |>
      mutate(race = case_when(
            race_1 == 1 ~ "White",
            race_2 == 1 ~ "Black",
            race_3 == 1 ~ "Hispanic",
            race_4 == 1 ~ "Asian",
            race_5 == 1 ~ "Other",
            race_5 == 1 ~ "Other",
            race_6 == 1 ~ "Other",
            race_8 == 1 ~ "Other"
      ),
      female = female_or_nb,
      interest = newsint, 
      pew_religimp = religiosity,
      pid3new = pid3,
      educ = case_when(
            educ == 1 ~ "No HS",
            educ == 2 ~ "HS degree",
            educ == 3 ~ "Some college",
            educ == 4 ~ "2-year",
            educ == 5 ~ "4-year",
            educ == 6 ~ "Post-grad",
            educ == 7 ~ NA_character_
      )) |>
      tbl_summary(type = list(interest ~ "continuous",
                              pew_religimp ~ "continuous",
                              pid3new ~ "categorical"),
                  label = list(female ~ "Female", 
                               age ~ "Age",
                               rural ~ "Rural",
                               educ ~ "Education",
                               pid3new ~ "Party ID",
                               race ~ "Race",
                               pew_religimp ~ "Religiosity (1-4)",
                               interest ~ "Political interest (1-4)"),
                  statistic = list(
                        all_continuous() ~ "{mean} ({sd})",
                        all_categorical() ~ "{p}%"),
                  missing = "no",
                  include = -c(starts_with("race_"), religiosity, 
                               newsint, pid3, female_or_nb)) |>
      modify_header(label = "Variable") |>
      bold_labels()
tbl2

##### S3 CHARACTERISTICS #####
tbl3 <- dcyber |> select(female, age, rural, educ, race, pid3new, interest, pew_religimp) |>
      mutate(pid3new = case_when(
            pid3new == 1 ~ "Democratic", 
            pid3new == 2 ~ "Independent",
            pid3new == 3 ~ "Republican",
      ),
      race = case_when(
            race == 1 ~ "White",
            race == 2 ~ "Black",
            race == 3 ~ "Hispanic",
            race == 4 ~ "Asian",
            race > 4 ~ "Other"
      ),
      educ = case_when(
            educ == 1 ~ "No HS",
            educ == 2 ~ "HS degree",
            educ == 3 ~ "Some college",
            educ == 4 ~ "2-year",
            educ == 5 ~ "4-year",
            educ == 6 ~ "Post-grad")) |>
      tbl_summary(type = list(interest ~ "continuous",
                              pew_religimp ~ "continuous",
                              pid3new ~ "categorical"),
                  label = list(female ~ "Female", 
                               age ~ "Age",
                               rural ~ "Rural",
                               educ ~ "Education",
                               pid3new ~ "Party ID",
                               race ~ "Race",
                               interest ~ "Political interest (1-4)",
                               pew_religimp ~ "Religiosity (1-4)"),
                  statistic = list(
                        all_continuous() ~ "{mean} ({sd})",
                        all_categorical() ~ "{p}%"),
                  missing = "no") |>
      modify_header(label = "Variable") |>
      bold_labels()
tbl3

# merged table
tbls <- tbl_merge(
      tbls = list(tbl1, tbl2, tbl3), 
      tab_spanner = c("**2021 CES**", "**2022 Prolific**", "**2023 CES**") 
)

as_kable(tbls, format = "latex")

#### APPENDIX FIGURES ####

##### PPS Distributions #####
# Study 1
drural |> mutate(rural.victims = ifelse(rural.victims == 1, 
                                        "Rural Victims", "Urban Victims")) |>
ggplot(aes(harm.economic, fill = factor(rural))) + 
      geom_histogram(size = 1, position = position_dodge(.1)) + 
      labs(x = "Perceived Problem Seriousness", 
           y = "Count") + 
      theme_bw() + 
      scale_fill_manual(values = c("darkgray", "darkgreen"), 
                         name = "R Identity",
                         labels = c("Non-rural", "Rural")) +
      facet_grid(~factor(rural.victims))

# Study 2
dpid |> filter(!is.na(condition_inparty_w_NP)) |>
      ggplot(aes(pps, color = factor(pid3))) + 
      geom_density(size = 1) + 
      labs(x = "Perceived Problem Seriousness", 
           y = "Density") + 
      scale_color_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      theme_bw() +
      facet_grid(~factor(condition_inparty_w_NP))

##### Policy Support Distributions #####

# Study 1
drural |> mutate(rural.victims = ifelse(rural.victims == 1, 
                                        "Rural Victims", "Urban Victims")) |>
      ggplot(aes(govt.spend, fill = factor(rural))) + 
      geom_histogram(size = 1, position = position_dodge(.1)) + 
      labs(y = "Count", 
           x = "Policy Support") + 
      theme_bw() + 
      scale_fill_manual(values = c("darkgray", "darkgreen"), 
                        name = "R Identity",
                        labels = c("Non-rural", "Rural")) +
      facet_grid(~factor(rural.victims))

# Study 2
dpid |> filter(!is.na(condition_inparty_w_NP)) |>
      ggplot(aes(policy.filters, fill = factor(pid3))) + 
      geom_histogram(size = 1, position = position_dodge(.1)) + 
      labs(x = "Policy Support", 
           y = "Count") + 
      scale_fill_manual(values = c("navyblue", "darkred"), 
                         name = "R Identity",
                         labels = c("Dem.", "Rep.")) +
      theme_bw() +
      facet_grid(~factor(condition_inparty_w_NP))


